In [1]:
import numpy as np
import pandas as pd
import ipywidgets as ipw
from ipywidgets import interact
from bokeh.models import ColumnDataSource
from bokeh.io import push_notebook, output_notebook, show
from bokeh.layouts import row, column
from bokeh.plotting import figure
output_notebook()
from collections import OrderedDict
old_settings = np.seterr(over='ignore') #Ignore warnings about overflow data points
This notebook serves to show how the Fermi-Dirac distribution and its derivative change with temperature.
The Fermi-dirac distribution is $$f(E)=\frac{1}{e^{(E-\mu)/k_bT}+1}$$
For the purpose of simplicity, we will plot the graph of $f(E-\mu)$. This allows us to avoid complications connected with $\mu$ being a function of temperature (consult the other notebook for this week or section 2.2.3 in the notes for more details)
We first define the relevant variables, such as temperature, $k_b$ and energy range. Then we calculate the distribution for each point on the graph and its derivative.
In [2]:
# Set up constant
N = 200; #number of point on graph
kB = 8.617e-5; #eV/K
T = 0.001; # Boltzmann constant measured interms of chemical potential
currenT = [f'T = {T:.2f} K']*N; #Create an array of labels to initiate interactive legend
kBT = kB * T;
E_min = -10.0;
E_max = 10.0; # Energy in eV
factor = 0.2;
# Set up the variable space for 1D plot
En_mu = np.linspace(E_min, E_max, N) #Energy measured in units of chemical potential
f = 1.0/(np.exp( (factor*En_mu)/kBT )+1) # Occupation number according to F-D statistics
df = np.gradient(f) # Take the first derivative of the F-D distribution.
The rest of the notebook deals with setting up the plots for the graphs
In [3]:
# Set up the first 1D plot (0K)
plot = figure(height=400, width=400, title="Fermi-Dirac statistics in 1D",
tools="pan,reset,save,wheel_zoom",
x_range=[E_min, E_max], y_range=[0.0, 1.1])
plot.xaxis[0].axis_label='Energy relative to chemical potential, E - \u03BC (eV)';
plot.yaxis[0].axis_label='Fermi-Dirac distribution f(E)';
source = ColumnDataSource(dict(x=En_mu, y=f, label=currenT));
r = plot.line(x='x',y='y',legend_group='label', source=source, line_color = 'blue', line_width=3, line_alpha=0.6);
In [4]:
# Set up the second interactive 1D plot
plot1 = figure(height=400, width=400, title="Fermi-Dirac statistics in 1D",
tools="pan,reset,save,wheel_zoom",
x_range=[E_min, E_max], y_range=[0.0, 1.1])
plot1.xaxis[0].axis_label='Energy relative to chemical potential, E - \u03BC (eV)';
plot1.yaxis[0].axis_label='Fermi-Dirac distribution f(E)';
source1 = ColumnDataSource(dict(x=En_mu, y=f, label=currenT));
rrrrr = plot.line(x='x',y='y',legend_field='label', source=source1, line_color = 'red', line_width=3, line_alpha=0.6);
In [5]:
# Set up a 1D plot of the first derivative of the F-D distribution at 0K
plot3 = figure(height=400, width=400, title="First derivative",
tools="pan,reset,save,wheel_zoom",
x_range=[E_min, E_max], y_range=[0.0, 0.6])
plot3.xaxis[0].axis_label='Energy relative to chemical potential, E - \u03BC (eV)';
plot3.yaxis[0].axis_label='First derivative df(E)/dE [1/eV]';
source3 = ColumnDataSource(dict(x=En_mu, y=-df, label=currenT));
rrrrrr = plot3.line(x='x', y='y', legend_group = "label", source=source3, line_color = 'blue', line_width=3, line_alpha=0.6);
In [6]:
# Set up an interactive 1D plot of the first derivative of the F-D distribution
plot4 = figure(height=400, width=400, title="First derivative",
tools="pan,reset,save,wheel_zoom",
x_range=[E_min, E_max], y_range=[0.0, 0.6]);
plot4.xaxis[0].axis_label='Energy relative to chemical potential, E - \u03BC (eV)';
plot4.yaxis[0].axis_label='First derivative df(E)/dE [1/eV]';
source4 = ColumnDataSource(dict(x=En_mu, y=-df, label=currenT));
rrr = plot3.line(x='x', y='y', legend_field = 'label', source=source4, line_color = 'red', line_width=3, line_alpha=0.6);
In [7]:
# Set up callbacks to update the 1D live graph
def update_data(T):
# Generate the new curve
kBT = kB*T;
currenT = [f'T = {T:.2f} K']*N;
En_mu = np.linspace(E_min, E_max, N) #Energy measured in units of chemical potential
rrrrr.data_source.data['y'] = 1.0/(np.exp( (factor*En_mu)/kBT )+1);
rrrrr.data_source.data['label'] = currenT;
rrr.data_source.data['y'] = -np.gradient(rrrrr.data_source.data['y']);
rrr.data_source.data['label'] = currenT;
push_notebook();
In [8]:
show(row(plot,plot3),notebook_handle=True);
In [9]:
interact(update_data, T = ipw.FloatSlider(min=0.001,max=1000.001,step=10,value=300,description = 'Temp (K)'));
interactive(children=(FloatSlider(value=300.0, description='Temp (K)', max=1000.001, min=0.001, step=10.0), Ou…
In [ ]:
In [ ]: